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Abstract 


We have undertaken the theoretical modelling of an existing liquid hydrogen (LH2) and gas hy- 
drogen (GH2) mixer subsystem of the El Ground Test Facility at NASA John C. Stenrds Space 
Center. The El test facility carries out comprehensive ground-based testing and certification of 
various liquid rocket engines and their components. The mixer described in this work is responsible 
for comb ining hi gh pressure LH2 and GH2 to produce a hydrogen flow that meets certain thermo- 
dynamic properties before it is fed into a test article. The desired properties are maintained by 
precise control of the mixture of LH2 and GH2 flows. The mixer is modelled as a general multi-flow 
l um ped volume for single constituent fluids using density and internal energy as states. The set of 
nonlinear differential equations is linearized about an equilibrium point and the resulting two-state, 
3-input linear model is analyzed as a possible candidate for control design. 


1 Introduction 


The Test and Engineering Directorate at NASA John C. Stearns Space Center (SSC) is continu- 
ing its efforts to assemble a software simulation environment that captures the static and dynamic 
characteristics of modem and future rocket engine systems. To that end, the design of a simulation 
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ing dy nami c fluid flow simulation and advanced control design. EASY/ROCETS is the fusion of 
EASY5x, a commercial package by Boeing [1], andROCETS (ROCket Engine Transient Simulator), 
a rocket engine analysis package developed by United Technologies uorporation under contract to 
NASA [2]. The user-friendly package is foreseen to fulfill the need for an accurate and verifiable 
rocket engine test facility simulation environment. 

We build upon the reports [3, 4, 5] and focus on one component of the high pressure hydrogen 
system of the El Test Facility, namely, the LH 2 and GH 2 mixer subsystem and associated flow 
controllers. The relevant components of the mixer subsystem are shown in the flow schematic of 
Figure 1. The LH 2 valve controls the flow of high pressure liquid J? 2 horn a pressurized tank; the 
GH 2 valve controls the flow of gaseous from a set of high pressure bottles; and the output valve 
controls the flow of if 2 into the test article. The primary objective of the mixer and control valves is 
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in spite of the steady depletion of the source bottles, measurement errors, modelling imperfections, 
and other uncertainties and perturbations. 


2 Background 

e Some application examples of the use of EASY/ROCETS are presented in [3] including FOR- 
TRAN listings of several modules and a brief user's manual. 

• The report [4] performs a detailed derivation of the model of the RUN-TANK module as 
well as its implementation within EASY/ROCETS. The low pressure and high pressure LOX 
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Figure 1: Diagram of HPH2 Mixer and Control Valves. 


and hydrogen subsystems are simulated with EASY/ROCETS and compared against previous 
computer simulations with good results. In addition, various control strategies were simulated 
for the high pressure hydrogen subsystem. The objective being the temperature regulation of 
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the lack of a suitable mathematical model of either the tank or the mixer. 


• In [5], several upgrades and enhancements to EASY/ROCETS were reported including the 
inclusion of the NIST-12 thermodynamic database. Also, a PID module is introduced based 
on the Allen-Bradley SLC 5/03 Programmable Logic Controller (PLC). During the course of 
this investigation, an inconsistency was discovered in the discrete-time version of the reported 
PID controller equations. The corrected equations will be presented in a separate report. 

• Most recently, Victor Marrero, one of the SSC 2001 Summer Undergraduate students, reported 
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hydrogen system with real test data. The results are very encouraging and prompts for 
continued analytical studies such as the one presented in this report. 


2.1 Results 

» We developed an analytical model of the mixer and control valves that is suitable for dynamic 
analysis and control design in the vicinity of a thermodynamic equilibrium point. 

* An optimal linear-quadratic regulator was designed and simulated using Matlab to illustrate 
the effectiveness of such a controller in compensating against perturbations that result in 
small deviations around the desired equilibrium. A simulation is included that indicates how 
one could study which valve or combination of valves is more effective in achieving a desired 
performance. In essence, this is a study on control authority or relegation. 
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s The linear model is compared with that obtained with Easy-5. This is a preliminary comparison 
but the results are very encouraging as many of the observed discrepancies are within tolerance. 


3 Mixer and Valve Models 

The mixer has been modelled within EASY/ROCETS with a two-input, one output VOLMOl 
block, which is a general multi-flow lumped volume for single constituent fluids using density and 
internal energy as states [3]. In this article, we will refer to the mixer’s internal energy and mass 
as the “state”, not to be confused with the definition of state of a system as described by its 
thermodynamic properties. Letting zi(f), z 2 (t ) denote the mixer’s specific internal energy and 
mass, respectively, the two dynamic equations are given by 

zi = ^ [U net - W net zj + Q] (1) 

Z2 = W net ( 2 ) 

where U ne t is the net (input minus output) energy into the mixer, W ne t is the net (input minus 
output) flow, and Q is the total heat transfer rate. The net flow is simply 

W net — w g + wi - w a 

where the subscripts stand for “gas”, “liquid”, and “out”, respectively. Since the block VOLMOl 
allows for flow reversals, then a total of eight possible expressions exist for the net energy U ne t 
depending on the direction (sign) of the input and output flows. In this study, only positive flows 
are considered, that is, the gas and liquid flows are into the mixer and the output flow is out of the 
mixer. In essence, normal operation requires that each valve always experience a positive pressure 
difference across it thus disallowing flow reversal. The net input energy term is then given by 

Unet = hgWg + hwi - h v w Q 

where h g , h u h v denote the enthalpy of the GH 2 stream, LH 2 stream, and volume block, respec- 
tively. 

The valve block within EASY/ROCETS is labelled VALV99 and it calculates the incompressible 
fluid flow value w f given the inlet P n y and outlet Pout pressures, the density of the stream p, and 
the valve flow coefficient Cff, as follows: 1 

Wf — \Jp p w {PiN ~ Pout ) Cff — apCff l 3 ) 

where p w is the density of water, and we have assumed that P/at > Pout ■ Substituting the 
expressions for the valve flows and the net internal energy and net flow into equations (l)-(2), and 

1 It is with very high probability that in order to achieve a closer match with physical reality, some of the valve blocks 
will have to be modified or replaced by a valve model that simulates compressible flow [6]. That study is left for future 
research. 
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treating the valve coefficients as inputs, the dynamic model of the mixer is a 2-state, 3-input system 
of nonlinear differential equations of the form 


z\ = Fi(zi, z 2 , Cf g , Cfi, Cf 0 ) (4) 

z 2 — F 2 {z\, z 2 , Cfg , Cfi, Cf 0 ) (5) 


where F\{-) and F 2 (-) are nonlinear functions of the state and valve coefficients. 

For the remainder of the article, we will denote constant or equilibrium values of any variable by 
an upper bar (•). Given constant values of valve flow coefficients Cf = [ Cf g Cfi Cf 0 ] T (superscript 
T denotes transposition) and constant fluid properties, the state of the model 


z{t) = 


z i : Internal Energy 

;?2 : Mass 


reaches a constant equilibrium point z. Next, consider perturbing such an equilibrium by small 
signals x(t) and u(t) so that 


z(t) = z + x(t) and Cf(t) = Cf + u(t ) 

where u(t) denotes a small valve coefficient correction signal. Then, a standard linearization of 
equations (4)- (5) results in the small signal model 

x — Ax + Bu (6j 


where x(t) is the small perturbation state vector, u(t) is the small perturbation control signal, and 
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where the partial derivatives are evaluated at the equilibrium state z and constant valve flow 
coefficient vector Cf. 
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under co nsi deration. Ideally, both x(t) and u(t) are zero; in practice, the state z(t) deviates from 
the desired z thus necessitating a corrective action by u(t). The model is suitable for a variety of 
analysis studies such as stability, controllability, and regulation/tracking of the equilibrium. The 
stability study is relevant since it reveals an important characteristic of the equilibrium z under no 
control, that is, with u(t) = 0. Appropriate action must be taken by the control u(t) in the event 
that z be unstable (worst case) or even “poorly behaved”. The notion of controllability reveals 
whether a suitable control u(t) exists such that the perturbation state x(t) may be steered to zero. 
Moreover, it is possible to study -which valve or combination of valves is most effective. Finally, the 
eq uili bri um regulation/tracking study refers to the actual control design to regulate x(t) to zero or 
to track a desired equilibrium trajectory as a function of time. 

Typically, in the control literature an output equation of the form 


y — Cx + Du 
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where matrices C and D are appropriately dimensioned is appended to the model (6) to account 
for the measurement of certain variables for example, temperature, pressure, or flow. The state 
postulate 2 implies that the state z(t) — z + x is directly available for control design if at least two 
independent thermo dynamic properties were measured, for example, temperature and pressure. 
Therefore, it is reasonable to assume at least initially that 

y = x = z — z (7) 


3.1 Calculation of Equilibrium Points 

By definition, an equilibrium point z of equations (4)-(5) satisfies 

zi = 0 =+> hgWg + hiidi - fhu ) a + <3 = 0 ( 8 ) 

z% — 0 =+ w g + wi - w 0 = 0 (9) 


Given a set of required output mass flow, pressure, and temperature w 0 , P 0 > T 0 , then the property 
routine MIPROPS was used to determine the enthalpy h a = h v . 3 Next, select an appropriate set 
of GH 2 and LH 2 properties. Then, using the equilibrium conditions (8)-(9) with <3 = 0, and the 
flow expressions (3) for each valve, the following two linear equations in the gas and liquid flow 
coefficients are obtained 
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(10) 


ini 


The mixer equilibrium state Z 2 is simply 


^2 = Vpv 

and z\ is determined from a suitable thermodynamic property table or from MIPROPS. 

2 For a simple compressible substance, any thermodynamic property is at most a function of two other independent 
thermodynamic properties. 

3 TP 


‘Enthalpy across each valve remains constant by conservation of energy. 
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To determine the linear perturbation model (6), it is necessary to find the analytic expressions of 
the indicated partial derivatives. For simplicity we include only those terms that do not evaluate 
to zero at the equilibrium point. Omitting the details, these are found to be: 
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Note that the evaluation of ^ and ^ would have to be done numerically within the EASY/ROCETS 
environment based on its property routines. Alternatively, we used the property software package 
called MIPROPS that asks the user to enter values for pressure, density, and temperature, and returns 
all the properties. A nominal property set corresponding to a mixer pressure of P v = 6000 psia and 
density p v = 2.914 lb/ ft 3 was chosen. By varying the pressure while keeping the density constant, 
we generated property data as shown in Table 1. The top plots of Figure 2 and the indicated 
straight-line fits provide the relation between enthalpy hv and pressure P v to the internal energy 
zi while keeping the mass z 2 (or equivalently the density) constant. The temperature that the 
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program converges to is seen to deviate slightly from the entered nominal value of T — 200 F: 

these deviations are likely to be caused by numerical convergence within the program itself. 

Similarly, data relating enthalpy h v and pressure P v to the mass £ 2 while keeping the internal 
energy z\ constant was obtained by entering a density value at a nominal temperature and adjusting 
the pressure until the returned value for the internal energy was z x = 407.4 BTU/lbm. This was 
repeated for a range of density values. Representative data is listed in Table 2 and the linear-fit 
relations are given in the bottom plot if Figure 2. 

Table 3 lists the pertinent equilibrium data corresponding to a given set of outflow requirements 
for a test-article. A second set of data is given in Table 4 corresponding to a new equilibrium that 
results when the properties of the gas and liquid in-streams change slightly. Note that the valve 
coefficients are not modified; consequently, the outflow properties deviate from the desired values. 

Table 5 lists nine data points of ( P , T, h) as a function of the independent variables Oi, z 2 ) for 
the fixed gas-, liquid-, and out-stream properties, and fixed valve positions given in Table 3. The 
equilibrium conditions (8) and (9) are computed and plotted in Figure 3 illustrating the fact that 
for a given set of in- and out-stream properties and fixed valve positions, the resulting equilibrium 
point is isolated. That is, no other equilibrium points exist within an infinitesimal neighborhood 
around that equilibrium point. We infer that is the case for any other equilibrium point. 


4 Control Design 


A preliminary control design has been completed and it is detailed in this section. Consider the 
linear model (6) valid in the neighborhood of the equilibrium listed in Table 3. The matrices A and 
B are evaluated to be 

609.65 -64.95 -33.99 

2.6001 2.1211 -0.6488 _ 

In classical control design one could proceed to work with the transfer functions associated with 
this system. Consider the two states as outputs, and since there are three inputs, then there are 
six transfer functions given by 


Ai(s) 

609.650 + 122.8) 

Ai(s) 

-64.950 + 406.6) 

U g (s) 

(s 4- 189.6)0 + 11.61) ! 

Ui(s) 

(s + 189.6)0 + 11.61) 

Xi(s) 

—33.990 + 8.86) 

*2 0) 

2.6(s - 104 ) 

Uo(s ) 

0 + 189.6) 0 + H.61) ’ 

U 9 (s ) _ 

0 + 189.6)0 + 11-61) 

*2(3) 

2.120 + 65.26) 

Ms) 

-0.65 0 + 12-25) 

Ui(s ) 

0 + 189.6)0 + 11.61) ’ 

Uo(s ) 

0 + 189.6)0 + 11-61) 


It is in general difficult to design PID-style controllers for each loop because of the interactions 
between channels. The approach chosen in this work is based on the state-space model (6) rather 
than on the individual transfer functions. 


-45.709 -7686.3 

-0.53844 -155.54 
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Given the system (6). a performance criterion of the form 
J = J°° (x T (t}Qx(t) + u T (t)Ru(t )) dt 

is minimized over all possible control signals u(t). For technical reasons, the matrix Q must be 
positive-semidefinite and penalizes how far the state x(t) deviates from the origin. Similarly, the 
matrix R must be positive-definite and penalizes how much control effort is used to drive the state 
x(t) to the origin. The resulting optimal controller is given by the linear full-state feedback law 

u(t) = — Kx(t) (12) 

where the matrix K is constant and is found by solving the steady-state Riccati Equation 

a t p + pa- pbr- 1 b t p + q = o =>• k = r- 1 b t p 

The closed-loop system 
x = [A — BK] x 

is guaranteed to be asymptotically stable, that is, x(t) 


■ 0 , t 


■ oo. 




' l 
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0 ' 

' 1 o' 

; R = 
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l 

0 

0 1 


0 

0 
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As a simulation example, select the following weighting matrices: 

Q = 

so that it is equally important to steer a* or x 2 to zero, and the three valves have equal authority 
in achieving the stated goal. The optimal gain matrix is found to be 

0.90578 -9.35771 

K = —0.14018 2.613 

-0.041332 0.182310 

The resulting control policy and state trajectories are shown in Figure 4 and Figure 5, respectively. 
Note that the gas valve shows to have the biggest influence. As a comparison, change tne control 

Weighting matrix R to 
20 0 0 


R 


0 

n 


1 0 
0 1 


so that the gas valve controller is penalized twenty times more than the liquid or out valves. The 
optimal gain matrix is found to be 
0.12697 -2.8667 

K = -0.54362 17.278 

-0.084231 0.85055 _ 

and the resulting control and state histories are given in in Figure 6 and Figure 7, respectively. 
Observe that indeed the authority of the gas valve has been reduced with a corresponding increase 
in the efforts of the liquid and out valves. Compared to the previous case, the perturbation states 
are seen to take slightly longer to reach the origin. 
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5 Conclusions and Further Research 


An analytical small-signal model of the high-pressure H 2 mixer and control valves has been derived. 
The model is suitable for a variety of analysis and control designs valid in the vicinity of an oper- 
ating or equilibrium point. A classical approach to control design would involve the consideration 
of a total of six single-input, single-output transfer functions. Unless the control loops are weakly 
coupled, such a design is difficult to complete because of the interaction among the control chan- 
nels. Instead, a state-variable control design has been pursued that results in an optimal full-state 
constant feedback law that can be implemented in real-time. As illustrated in the simulations, it 
is straightforward to design the regulator and it is possible to assign different authorities to each 
control valve. 

Topics for further research include: 

• Introduce a valve model for compressible flows as well as valve dynamics, such as friction and 
possibly other nonlinear effects. It has been experimentally found that there is a “dip” in the 
“Flow” versus “Valve Position” curve that may result in unstable operation [7]. A heuristic 
approach has been used to solve this problem. 

• Investigate the problem of equilibrium trajectory tracking. This refers to the design of control 
policies that steer the PVT - curve as a function of time. 


Investigate the combined tank/mixer control problem. This is a six-state and four control- valve 
problem. Clearly, the quality of the out-flow properties is affected by the precise regulation of . 
both the tank and mixer internal dynamic states. 

Implement these findings within EASY/ROCETS. The following represents a modest effort 
along these lines. The isolated mixer has been built in EASY/ROCETS as shown in Figure 8. 
The “Linear Analysis” module was used to identify the resulting equilibrium point yielding 
p v = 1 .686481?— 03 Ibm/in 3 and E v = 407.32 BTU/lbm. These compare very well with z 2 /V = 
1.685E - 03 Ibm/in 3 and = 407.4 BTU/lbm, respectively. The same module generates the 
following linearized system model (the subscript ER denotes the EASY/ROCETS matrices): 


Aer 


-60.34 -7634.57 


-0.919 


-155.7 


Ber 


608.69 -64.82 
2.59 2.12 


-33.97 

-0.649 


which are to be compared with the analytical matrices 

[ -45.709 -7686.3 

[ -0.63844 -155.54 

The cause for the discrepancy in the first column of matrices Aer and A should be investigated 
further. The rest of the entries are deemed to be within acceptable tolerances. Finally, the 
system’s natural frequencies are the eigenvalues of the system matrix Aer ■ These are easily 
found to be -11.39, - 203.7. The eigenvalues of A are -11.61, - 189.6. The natural or 
unforced responses corresponding to these sets of eigenvalues (poles) are quite similar. 

• Last but definitely not least , we must validate these findings with real data. 


609.65 -64.95 -33.99 

2.6001 2.1211 -0.6488 
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P (psia) 

T(F) 

zl (BTU/lb) 

h (BTU/lb) 

5980 

-200.7 

405.2 

785.2 

5985 

-200.51 

405.8 

786.1 

5990 

-200.33 

406.3 

787 

5995 

-200.15 

406.9 

787.8 

6000 

-199.96 

407.4 

788.7 

6005 

-199.78 

408 

789.6 

6010 

-199.6 

408.5 

790.5 

6015 

-199.41 

409.1 

791.3 

6020 

-199.23 

409.6 

792.2 1 


Table 1: Hydrogen Properties at Constant Density p — 2.914 lb/ ft 3 . 


P (psia) 

T(F) 

h (BTU/lb) 

P W, ft 3 ) 

5595 

-200.88 

777.4 

2.8 

5770 

-200.47 

782.3 

2.85 

5949 

-200.08 

787.3 

2.9 

5985 

-200.01 

788.2 

2.91 

6000 

-199.96 

788.7 

2.914 

6132 

-199.7 

792.3 

2.95 

6320 

-199.32 

797.5 

3 

6513 

-198.94 

802.8 

3.05 


Table 2: Hydrogen Properties at Constant Internal Energy z\ = 407.4 BTU/lb. 
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h 

fqn-’Tr/ihrp) 

C'vopen 

%OPEN 

w 

77vm I q 

Z\ 

BTU/ibm 

^2 

Ibm 

gh 2 in 

13500 

90 

2.91 

2113.6 

230 

4.35 

- 

- 

- 

gh 2 OUT 

6000 

135 

1.53 

2113.6 

- 

- 

26.0 

- 

gap 

lh 2 in 

8500 

-400 

5.81 

184.6 

115 

23.4 

- 

- 

B 

lh 2 OUT 

6000 

-370 

5.02 

184.6 

- 

- 

57.0 

- 

1 

Mixer 

6000 

-200 

2.91 

788.7 

- 


- 

407.4 

7.28 

Outflow 

5533 

-198 

2.76 

788.7 

270 

47.37 

83.0 

- 

- 


Table 3: Equilibrium Data. Underlined Values are Computed. (IN=Into Valve; OUT=Out of Valve) 
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Figure 2: Top: if 2 Enthalpy and Pressure Variation with Internal Energy at Constant Density. Bottom 
i ?2 Enthalpy and Pressure Variation with Mass at Constant Internal Energy. 
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Table 4: Perturbed Equilibrium Data. Underlined Values are Computed. (IN— Into Valve; OUT— Out 
of Valve) 
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Table 5: Data Used to Compute the Equilibrium Conditions EQl (8) and EQ2 (9). 
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Gas Valve Small Control 
(Authority 1) 



Figure 4: Optimal Full-State Feedback Controllers for Case 1 



Figure 5: Optimal State Trajectories for Case 1. 
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Liquid Valve Small Control 
(Authority 1) 



Figure 6: Optimal Full-State Feedback Controllers for Case 2. 



Figure 7: Optimal State Trajectories for Case 2. 
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VALV99 from ROCETS- 



Mixer(10) 


Figure 8: EASY/ROCETS Implementation of the Mixer and Control Valves. 
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